Dynamical mean-field study of the Mott transition in thin films 
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The correlation-driven transition from a paramagnetic metal to a paramagnetic Mott-Hubbard 
insulator is studied within the half-hlled Hubbard model for a thin-film geometry. We consider 
simple-cubic films with different low-index surfaces and film thickness d ranging from d = 1 (two- 
dimensional) up to d = 8. Using the dynamical mean-field theory, the lattice (film) problem is self- 
consistently mapped onto a set of d single-impurity Anderson models which are indirectly coupled 
via the respective baths of conduction electrons. The impurity models are solved at zero temperature 
using the exact-diagonalization algorithm. We investigate the layer and thickness dependence of the 
electronic structure in the low-energy regime. Effects due to the finite film thickness are found to be 
the more pronounced the lower is the film-surface coordination number. For the comparatively open 
sc(lll) geometry we find a strong layer dependence of the quasi-particle weight while it is much less 
pronounced for the sc(110) and the sc(100) film geometries. For a given geometry and thickness d 
there is a unique critical interaction strength U C 2(d) at which all effective masses diverge and there 
is a unique strength U c i(d) where the insulating solution disappears. U C 2(d) and t/ c i(d) gradually 
increase with increasing thickness eventually approaching their bulk values. A simple analytical 
argument explains the complete geometry and thickness dependence of U C 2- U c i is found to scale 
linearly with U C 2- 



I. INTRODUCTION 

Electron-correlation effects in systems with thin-film 
geometry have gained increasing interest in condensed- 
matter physics. In particular, there has been in- 
tense research on thermodynamic phase transitions to a 
symmetry-broken (e. g. magnetic) state below a critical 
temperature jl], g{ . In magnetic thin films the thickness 
dependence of the order parameter, of the critical tem- 
perature as well as of the critical exponents has been 
investigated both, experimentally J3|, |^, ||, || and theo- 
retically g | g @. 

A thin film in three dimensions belongs to a two- 
dimensional universality class, regardless of the film 
thickness d Due to the Mermin- Wagner theo- 

rem [fl2|| , however, an effectively two-dimensional spin- 
isotropic system cannot display long-range magnetic or- 
der at any finite temperature. This is one important 
reason why anisotropies play a fundamental role for the 
understanding of thermodynamic phase transitions in 
thin films. The necessary inclusion of anisotropies, how- 
ever, makes a theoretical description considerably more 
complicated. 

For a quantum phase transition the situation is differ- 
ent: Symmetry breaking need not occur at the transition 
point, and the energy scale that is characteristic for the 



transition at zero temperature, remains meaningful at 
any finite temperature. Consequently, anisotropies are 
not vital for the understanding of a quantum phase tran- 
sition in thin films: The transition can be studied within 
an isotropic model and at any temperature, starting 
from the monolayer (d = 1) up to the three-dimensional 
limit (d <— > oo). From this point of view, following up 
the characteristics of a quantum phase transition as a 
function of d, may be the better defined and the simpler 
problem if compared with a thermodynamic transition. 

One of the prime examples for a quantum phase tran- 
sition is the correlation-driven transition from a para- 
magnetic metal to a paramagnetic insulator jl3|, [l4|. 
Generally, the Mott transition is of interest since strong 
electron correlations lead to low-energy electronic prop- 
erties that cannot be understood within an independent- 
electron picture; the important correlation effects must 
be treated non-perturbatively. Conventional band the- 
ory is unable to provide a satisfactory description of the 
transition. 

While magnetic phase transitions in thin films are 
traditionally studied within localized-moment models 
(Ising or Heisenberg model) [Q, [§, ^, |lO), the possibly 
simplest generic model for the Mott transition is the 
Hubbard model jl5| . Contrary to the localized-moment 
models, the Hubbard model describes itinerant electrons 
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on a lattice which may form local moments as a conse- 
quence of the strong on-site Coulomb interaction. Start- 
ing from the early approaches of Mott Hubbard 
fHI , and Brinkman and Rice JlTj , there has been exten- 
sive work on the Mott transition in the Hubbard model 
(for a recent overview see Ref. [fu]). A thin-film geom- 
etry has not been considered up to now. 

The following particular questions shall be addressed 
in the present study of the thin-film Mott transition: 
First, the breakdown of translational symmetry in the 
film normal direction introduces a layer dependence of 
physical quantities that characterize the transition. The 
layer dependence of the effective mass m* or the dou- 
ble occupancy (n-|-n^), for example, is worth studying. 
Second, it has to be expected that there is a depen- 
dence on the film geometry. We can distinguish between 
surface effects and effects due to the finite film thick- 
ness: Surface effects are already present in thick films 
(d i— » oo) and will be more pronounced for films with 
comparatively "open" surfaces, i. e. for surfaces where 
the nearest-neighbor coordination number is strongly re- 
duced. For thin films there may be finite-thickness ef- 
fects in addition. Here, the perturbation of the elec- 
tronic structure introduced by one of the film surfaces 
affects the electronic structure in the vicinity of the 
other one; both surfaces can "interact" with each other. 
Third, it is not clear from the beginning whether or not 
there is a unique critical interaction strength U c at which 
the whole film undergoes the transition. Analogously to 
magnetic transitions where an enhanced surface critical 
temperature is discussed [l9[ ^0|, a modified U c for 
the film surface might exist for thicker films. Finally, 
it is interesting to see the critical interaction strength 
U c evolving as a function of d and to investigate how it 
approaches the bulk value. Again, the crossover from 
D = 2 to D = 3 dimensions should depend on the par- 
ticular film geometry. 

In the recent years, comprehensive investigations of 
the Mott transition have been performed (l4|, pil|] fo r the 
Hubbard model in infinite spatial dimensions |22|, p3| . 
The D = oo model is amenable to an exact solution 
by a self-consistent mapping onto an effective impurity 
problem which, however, must be treated numerically. 
Neglecting the spatial correlations, the same method can 
be applied for an approximate solution of the Hubbard 
model in any finite dimension D < oo. This constitutes 
the so-called dynamical mean-field theory (DMFT) |25| , 
|2l| which will be employed also for the present study. 
It is another intention of the paper to demonstrate that 
DMFT can successfully be applied to a film geometry. 

The mean-field treatment of the Mott transition in 



D < oo rests on the assumption that spin- and charge 
fluctuations are reasonably local. In particular, the elec- 
tronic self-energy is a local quantity within DMFT , 
(E) = 5ijT,i(E). The relevance of non-local contribu- 
tions for a qualitatively correct description of the Mott 
transition is not well understood at present. Effects of 
the non-locality of the self-energy can be estimated to 
some extent by conventional second-order U perturba- 
tion theory (SOPT). As is shown in Ref. [p5| , the local 
approximation is rather well justified for a D — 3 simple- 
cubic lattice. Non-local contributions become more im- 
portant for D = 2 and especially for D = 1. With 
respect to the film geometry, the question is whether or 
not non-local contributions can be neglected also near 
the film surfaces. This has been investigated recently 
by means of SOPT applied to semi-infinite simple-cubic 
lattices j26|. The result is that at the surface the local 
approximation is as well justified as in the bulk. 

There is an additional suppression of the effects of 
non-local fluctuations at half-filling: The low-energy 
electronic structure, being relevant for the transition, 
is governed by the (real) linear expansion coefficient 
fiij = dT,ij(E) j 'dE\E=o- Since ReEij(-E) is a symmet- 
ric function of E in the particle-hole symmetric case 
and for nearest neighbors i and j, we have = 0. 
Within SOPT the second nearest-neighbor term fyj can 
be shown |2(| to be smaller by about two orders of mag- 
nitude compared with the local one. 

This shows that - at least for weak coupling - a local 
self-energy functional is a reliable approximation. Be- 
yond the weak-coupling regime, however, the approxi- 
mate locality actually is an assumption whose appro- 
priateness has not yet been verified. We nevertheless 
expect DMFT to be a good starting point to study the 
Mott transition in a D = 3 film geometry. 

We start our investigations by specifying the film 
geometries to be considered and discuss the DMFT 
approach with respect to thin films. In the first part 
of Sec. Ill we briefly focus on the Mott transition in 
the infinitely extended and translationally invariant 
model. After that the results for the films are analyzed 
in detail. Sec. IV concludes the study. 

II. DMFT FOR HUBBARD FILMS 

We consider the Hubbard model at half-filling and 
zero temperature. Using standard notations, the Hamil- 
tonian reads: 

H = Yl *i> C L<> + n il7 ni- a . (1) 



2 



The hopping integrals £y are taken to be non-zero be- 
tween nearest neighbors i and j. t = —tuj) = 1 sets 
the energy scale. The thin-film geometry is realized by 
assuming i and j to run over the sites of a system that is 
built up from a finite number d of adjacent layers out of 
a three-dimensional periodic lattice. We study simple- 
cubic films with surface normals along the low-index 
[100], [110] and [111] directions. The model parameters 
are taken to be uniform, i. e. t and U are unchanged at 
the two film surfaces. With respect to the film nor- 
mal, translational symmetry is broken; lateral trans- 
lational symmetry, however, can be exploited by two- 
dimensional Fourier transformation: The hopping tij is 
transformed into a matrix e a( 3(k) which is diagonal in 
the wave vector k of the first two-dimensional Brillouin 
zone (2DBZ). a, (3 — 1, . . . , d label the different layers. 
For nearest-neighbor hopping and low-index sc films the 
Fourier-transformed hopping matrix is tridiagonal with 
respect to the layer indices. Its non-zero elements are 
given by: e||(k) = e QQ (k) and ej_(k) = e aa±1 (k). Let a 
denote the lattice constant. For a sc(100) film the lateral 
and normal dispersions then read: 

eii(k) = 2t(cos(k x a) + cos(k y a)) , 

e ± (k) 2 = t 2 . (2) 

Only the square of e^(k) will enter the physical quan- 
tities we are interested in. For the sc(110) geometry we 
have: 

eii (k) = 2t cos(k x a) , 

e_L(k) 2 = 2t 2 + 2t 2 cos(V2k y a) , (3) 
and the dispersions for the sc(lll) films are: 
c||(k) = 0, 

e_L(k) 2 = M 2 + 2t 2 cos(V2k y a) 

+ At 2 cos{yfJj2k x a) cos{^/lj2k y a) .(4) 

Note that the parallel dispersion vanishes since there 
are no nearest neighbors within the same layer. Due 
to the finite film thickness d, the local free (U = 0) 
density of states acquires a layer dependence: pf 1 ^ (E) = 
Pq ^ (E) for sites i within layer a. Particle-hole symmetry 
requires pi (E) to be a symmetric function of the energy 
for each layer: On the bipartite lattice the odd moments 
^ l+l pS\E)dE (I = 0,1, .. .) vanish for all a (cf. Ref. 

To study the transition from a paramagnetic metal 
at weak coupling to a paramagnetic Mott-Hubbard in- 
sulator at strong U, we restrict ourselves to the spin- 
symmetric and (laterally) homogeneous solutions of the 



mean-field equations. As usual |14| we thereby ignore 
antiferromagnetic ordering which is expected to be re- 
alized in the true ground state at any U > 0. The on- 
site Green function Gu{E) = {{c ia \ cJ ct ))e thus depends 
on the layer index only: Ga(E) = G a (E). The same 
holds for the self-energy T*n(E) = T, a (E) which is a 
local quantity within the mean-field approach. Via two- 
dimensional Fourier transformation we obtain from the 
Dyson equation: 

II k 

i? a/ 3(k, E) = {E + n)5 a p - e Q/ 3(k) - 8 a p?, a {E) , 

(5) 

where N\\ is the number of sites within each layer (jVii i— > 
oo) and k e 2DBZ. For the particle-hole symmetric 
case, the Fermi energy is given by /i = U/2. Since e Q /3(k) 
is tridiagonal, the matrix inversion is readily performed 
numerically by evaluating a continued fraction of finite 
depth which is given by the film thickness d (cf. Ref. 
f§)- 

The layer-dependent self-energy £ Q (E) shall be calcu- 
lated within the dynamical mean-field theory (DMFT) 
|25| , pl[ . This non-perturbative approach treats the lo- 
cal spin and charge fluctuations exactly. Neglecting 
spatial correlations, a homogeneous lattice problem can 
be mapped onto an effective impurity problem supple- 
mented by a self-consistency condition j2^, |29| . For the 
present case of a thin Hubbard film we have to account 
for the non-equivalence of sites within different layers. 
Therefore, the film problem is mapped onto a set of d 
different impurity models with d self-consistency condi- 
tions. 

The DMFT equations are solved by means of the fol- 
lowing iterative procedure: We start with a guess for the 
self-energy T, a (E). This yields the on-site Green func- 
tion G a (E) via Eq. (||). In the next step we consider 
a single-impurity Anderson model (SIAM) |30| for each 
layer a: 

a <r,fc=2 

+ E (^U+H.c.) . (6) 

a,k=2 

Here ej^ and denote the conduction-band energies 
and the hybridization strengths of the a-th SIAM, re- 
spectively. It is sufficient to fix the free (U = 0) impurity 
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Green function G^ (E) which is obtained from the a-th 
DMFT self-consistency condition as: 



GW(E) = (G a (E)-i + X a (E)y 



(7) 



The crucial step is the solution of the impurity models 
for a = 1,. .. ,d to get the impurity self-energy T, a (E) 
which is required for the next cycle. 

The computational effort needed for the solution of 
the impurity models scales linearly with the system size. 
It is enhanced by a factor d/2 compared with DMFT 
applied to the translationally invariant (bulk) Hubbard 
model if one takes into account the mirror symmetry 
with respect to the central layer of the film. Com- 
pared with the translationally invariant problem, we 
only found a slight increase in the number of cycles nec- 
essary for the convergence of the iterative procedure. 
The coupling between the different impurity problems 
via their respective baths of conduction electrons turns 
out to be weak. 

The application of DMFT to the Hubbard model in 
thin-film geometry rests on exactly the same assump- 
tion that is necessary for the application of DMFT to 
any finite-dimensional system, namely on the local ap- 
proximation for the self-energy. To be precise: the self- 
energy is taken to be local, £^(£7) — 8ijHi(E), and to 
be given by the (diagrammatic) functional of the full 
local propagator Gu(E) only. This is sufficient to es- 
tablish the mapping onto the impurity models and to 
derive the self-consistency condition. Since the local ap- 
proximation is the only approximation used so far, the 
hopping between the layers is treated on the same level 
as the hopping within each layer. Near the film center 
and in the limit of infinite film thickness, our approach 
thus recovers the D = 3 (sc) bulk properties, irrespective 
of the particular film surface geometry. This represents 
a non-trivial check of the numerics. 

For the solution of the SIAM we employ the exact- 
diagonalization (ED) method of Caffarel and Krauth 
pljj . The ED has proven its usefulness in a number of 
previous applications [|T[ |33. 34, 2l[. The main idea 
is to consider a SIAM with a finite number of sites n s . 
This results in an obvious shortcoming of the method: 
ED is not able to yield a smooth density of states. An- 
other disadvantage consists in the fact that finite-size 
effects are non-negligible whenever there is a small en- 
ergy scale relevant for the problem to be investigated. 
With respect to the Mott transition, finite-size effects 
become important close to the critical interaction. On 
the other hand, there are a number of advantages: The 
ED method is based on a simple concept; it is easy 



to handle numerically and computationally fast if com- 
pared with the quantum Monte-Carlo (QMC) approach 
H|, m ||^, H|] . This is of crucial importance for a 
systematic study that covers a large parameter space. 
Opposed to QMC, the ED method is particularly well 
suited to study the model at T = 0. 

Within the ED algorithm the functional equations 
for the self-energy are solved on a discrete mesh on 
the imaginary energy axis: iE n = i(2n + 1)tt//3 (n = 
0, 1, . . . , n max ). The fictitious inverse temperature (3 
defines a low-energy cutoff, n ma x represents a high- 
energy cutoff. Eq. (Q) provides the bath Green function 
Ga\iE n ) from which we have to fix the parameters of a 



SIAM with n s sites. Following Ref. f31 
by minimization of the cost function 



this is achieved 
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1 J2\G { ° ) (*E n )- 1 -G£l(iE n 



71=0 



(8) 



(a) 



with respect to the conduction-band energies e fc 

and the hybridization strengths Vu (k = 2, ...,n s ). 

Thereby, Ga\iE n ) is approximated by the free (U = 0) 
Green function of an n s -site SIAM: 



G£) (iE n -fx)- 



iE„ 



E 

k=2 



0)\2 



iE n 



» 



(9) 



Obviously, the method is exact for n s i— > 00 only. The 
convergence with respect to n s , however, has been found 
to be exponentially fast fl3l| . Typically n s — 6 — 10 sites 
are sufficient for interaction strengths not too close to 
the critical interaction. The results are independent 
of the low-energy cutoff provided that is chosen 
to be sufficiently small. Errors show up in the critical 
region close to the transition. Compared with the error 
due to the finite n s , however, this is negligible. Once 
the SIAM is specified, Lanczos technique |27j may be 
employed to calculate the ground state and the T = 
impurity Green function G a (iE n ). The local self-energy 
of the a-th layer can be derived from the impurity 
Dyson equation T, a (iE n ) = Ga } (iS„) _1 - Gi mp («£'n) _1 . 



III. RESULTS AND DISCUSSION 



A. Bulk 

Let us first concentrate on the T = Mott tran- 
sition in the translationally invariant Hubbard model 
before we come to the discussion of the results for 
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Fig. 1. Results for the (bulk) Mott transition in the half- 
filled Hubbard model as obtained within DMFT-ED. U de- 
pendence of the linear coefficient f3 in the low-energy ex- 
pansion of £(.E) for the metallic (solid line) and the insu- 
lating solution (dashed line). The critical interactions U c i 
and U C 2 are indicated by the vertical lines. Lower panel: 
1/E coefficient a in the low-energy expansion for the insu- 
lating solution. Calculation for a D = 3 simple-cubic lattice. 
Nearest-neighbor hopping t = 1. Width of free density of 
states: W = 12. n s = 8. 



the film geometry. This case has been the subject 
of numerous DMFT studies during the recent years 
H, H H H ||, HI HI H |§ . Most investigations re- 
fer to the Bethe lattice with infinite connectivity where 
we have a semi-elliptical free density of states. However, 
within DMFT no qualitative changes are expected when 
considering the D — 3 simple-cubic lattice which also 
yields a symmetric and bounded free density of states. 
The D = 3 sc lattice is considered here since it repre- 
sents the limit of infinite film thickness d i— > oo with 
respect to our film results. 

Fig. 1 shows the U dependence of the coefficients a 
and (3 in the low-energy Laurent expansion of the self- 



/ , a U „ 
^(E) = - + -+(3E- 



(10) 



The calculation is performed for n s = 8. For a metal, 
i. e. if a = 0, the coefficient (3 yields the usual quasi- 
particle weight z = (1 — |^3|. We find a metallic 

solution for interaction strengths up to a critical value 
U c2 = 16.0 (U e2 = AW/3 in terms of the width of the 
free density of states W — 12). As U approaches U C 2, 
(3 diverges, i. e. the quasi-particle weight vanishes. At 
U C 2 the metallic solution continuously coalesces with the 
insulating solution that is found for strong U. For de- 
creasing U the insulating phase ceases to exist below 
another critical interaction strength U c i = 11.5 which is 
marked by the vanishing l/E expansion coefficient a as 
well as by (3 i— > — oo (Fig. 1). We find t/ cl < U c2 ; there 
is a region where both, the metallic and the insulating 
solution, coexist Q. 

A precise determination of the critical interactions, 
at least of <7 c2 , is not possible by means of the ED 
method (see next section). Qualitatively, however, the 
results are consistent with the findings for the D = oo 
Bethe lattice || |l|, || |J, H EH IS PL Within 
the iterative perturbation theory (IPT) fr28|, [2l[ a nar- 
row quasi-particle resonance is seen to develop at the 
Fermi energy for increasing U in the metallic solution. 
The spectrum has a three-peak structure, two additional 
charge-excitation peaks (Hubbard bands) show up at 
E ~ ±U/2. For U >—> U c2 the effective mass z~ x diverges 
as in the Brinkman-Rice variational approach jTjJ. On 
the other hand, at strong U the Hubbard bands are well 
separated by a gap in the insulating solution as in the 
Hubbard-III approach (T^|. One can follow up the in- 
sulating solution by decreasing U down to U = U c \. 
For T = a coexistence of the solutions (U c i < U C 2) 
has first been observed within IPT |}9| |4(J. It is con- 
firmed by ED |53f as well as by means of a recent numer- 



ical renormalization-group calculation (NRG) J42|, 45 
which is particularly suited to study the critical regime. 
The comparison between the respective internal ener- 
gies within IPT HI, ED HI |U and NRG El as well 



as a simple argument mentioned in Ref. [}41| show that 
the metallic phase is stable against the insulating one 
in the whole coexistence region up to [/ c2 : The T = 
transition is found to be of second order. 

One should also be aware of serious physical ar- 
guments [^6|, [l4| which have been raised against a 
transition scenario with two different critical interaction 
strengths. The exhaustion problem mentioned in Ref. 
ptj] at least shows that a conclusive understanding 
of the Mott transition in D = 00 has not yet been 
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achieved. From the above discussion and our own 
results, however, we can conclude that the numerical 
evidences for the existence of a finite coexistence region 
are strong. Another problem is tackled in Ref. Jt7| , 
where the concept of a preformed gap Q] is shown to 
be at variance with Fermi liquid theory. Our ED study 
cannot contribute to settle this interesting question 
since the detailed picture of the low-energy electronic 
structure in the limit U i— » U c i is concerned. 



B. Films 

In the following we discuss our results obtained by 
the ED approach to investigate the characteristics of the 
Mott transition in thin Hubbard films. Routinely, the 
calculations for the Hubbard films have been performed 
with n s = 8 sites in the effective impurity problems. We 
systematically compared with n s — 6 and also checked 
against n s — 10 at a few data points. It turns out that 
there are no significant differences between the results 
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Fig. 2. Imaginary part of the layer-dependent self-energy 
Eq. on the discrete mesh of the imaginary energies iE n = 
(2n + 1)tt/3 ((3 = 50, /T 1 = 0.0016 W) for the d = 5 sc(100) 
film at U = 14. Results for the three inequivalent layers 
a = 1 — 3 as indicated, a = 1: surface layer. Solid lines: 
metallic phase. Dashed lines: metastable insulating phase. 
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Fig. 3. The imaginary part of the Green functions corre- 
sponding to the self-energies in Fig. 2. 



for the different n s as long as U is not too close to U C 2 
(see Fig. 10 and the related discussion). Choosing a 
small fictitious temperature = 0.0016 W and a large 
high-energy cutoff (2n max — l)vr//3 > 2U ensures the re- 
sults to be independent of the discrete energy mesh (see 
also discussion of Fig. 10). For the electron-hole sym- 
metric case we can reduce the number of parameters in 
the multi-dimensional minimization [Eq. (^)] and for the 
DMFT self-consistency [Eq. (§)] by setting e k = -e k > 



and V? = V, 2 , for k, k' with k + k' = 



2. For n, = 



and d = 6, which are typical values considered here, 
there are (n s — l)d/2 = 21 parameters to be determined. 
The stabilization of paramagnetic solutions has turned 
out to be unproblematic. A mixing of "old" and "new" 
parameters (50%), however, has been found to be neces- 
sary to obtain a converging self-consistency cycle for the 
sc(lll) films. Apart from the coexistence of a metal- 
lic and an insulating solution in a certain U range, we 
always found a unique solution of the mean-field equa- 
tions. 

Fig. 2 shows the imaginary part of the self-energy 
as obtained for a d — 5 sc(100) film. At U = 14 we 
find two solutions. In the metallic one there is a signifi- 
cant layer dependence of Im S Q (iE) with a considerably 
larger slope dY^ a jdE for E ^ in the surface layer 
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(a = 1). This is plausible since the surface-layer sites 
have a reduced coordination number n resulting in a di- 
minished variance A = nt 2 of the free local density of 
states. Thus U/y/A is larger compared with the film 
center tending to enhance correlation effects. The large 
values of dH a jdE indicate that the system is close to 
the transition. 

At E = the imaginary part of the self-energy van- 
ishes for all layers as it must be for a Fermi liquid (note 
that in Fig. 2 Im X is plotted on the discrete mesh iE n 
only). A weak layer dependence is noticed for the in- 
sulating solution. In this case ImE Q (i£) diverges for 
E i ► . 

The imaginary parts of the corresponding Green func- 
tions are shown in Fig. 3. Again, the layer dependence is 
more pronounced for the metallic solution. For the insu- 
lating solution we have ImGo,(iO + ) = 0, i. e. a vanishing 
layer density of states at E = 0. Contrary, ImG Q (iO + ) 
stays finite for the metal. The E = value is given by 
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Fig. 4. Momentum distribution function n a (k) for the (non- 
equivalent) layers a = 1,2,3 along a high-symmetry direc- 
tion in the irreducible part of the 2DBZ. Results for the d — 5 
sc(100) film at different U. 



Eq. (||) where 

i? Q/3 (k,z0 + ) = 



(/i + i0 + )d a/3 - eapQt) - <5 Q/3 £ Q (0) 
i0 + 5 a p - e a p(k) = R a p(k, i0 + ) 



(7=0 

(11) 



and where we used that /i = S a (0) = (7/2 for all 
layers a. This implies that (for a local Fermi liq- 
uid) the layer density of states at the Fermi energy 
/Oq(0) = — Im G a (i0 + )/ it is unrenormalized by the inter- 
action. As for [7 = 0, however, it is layer dependent. For 
the translationally invariant Hubbard model with local 
self-energy the pinning of the density of states is a well 
known and general consequence of Luttinger's sum rule 
H). In the case of Hubbard films (and within DMFT) 
there is a pinning of p a (0) only for the case of half-filling 
since off half-filling £ Q (0) acquires a layer dependence. 

From the low-energy expansion of the self-energy for 
the metallic solution we can calculate the so-called layer- 
dependent quasi-particle weight 



dZ a (E) 



dE 



E=0 



(12) 



Once a self-consistent solution of the mean-field equa- 
tions on the discrete energy mesh iE n has been obtained, 
the self-energy can be determined in the entire complex 
energy plane, and there is no difficulty to calculate the 
derivative in Eq. (|l2|). Let us briefly discuss the physical 
meaning of z a which for a film geometry is slightly dif- 
ferent compared with the bulk case [Q. Exploiting the 
lateral translational symmetry and performing a two- 
dimensional Fourier transformation, the Green function 
at low energies is given by: 



y/zp. (13) 



1 II 3 II 



El - T(k) 



a/3 



ZnCapityy/zp is the renormalized 
in and ju label the N\\ sites in the 



Here (T(k)) Q)9 
hopping matrix 
layers a and /3, respectively, and k £ 2DBZ. Only 
the linear term in the expansion of the self-energy is 
taken into account. For each wave vector k the matrix 
T(k) can be diagonalized by a unitary transformation 
which is mediated by a matrix J7 ar (k). The d eigenval- 
ues T} r (k) of T(k) (r = 1, . . . , d) yield the dispersions of 
the quasi-particle bands near the Fermi energy and de- 
termine the d one-dimensional Fermi "surfaces" in the 



7 



2DBZ via ?7 r (k) = 0. As can be seen from Eq. (|13[), 
when k approaches the r-th Fermi surface there is a 
discontinuous drop of the a-th momentum-distribution 
function 

1 f° 

n a (k) = / ImG aa (k,E + iO + )dE , (14) 

^ J-oo 

which is given by: 

£n Q (k F ) = |[/ Qr (kp)| 2 • z a . (15) 
Summing over the d Fermi surfaces, we have: 

d 

z a = y^fa a (k F ) , (16) 

r=l 

The momentum-distribution function can be calcu- 
lated on the imaginary energy axis via: 



d=5 



1 2 °° 

(k) = - + lim -ReV G aa (k, iE n 



(17) 



Fig. 4 shows n Q (k) along a high-symmetry direction in 
the 2DBZ for the d = 5 sc(100) film. For the non- 
interacting system the 5 bands ?7r(k) are occupied at 
the r = (0,0) point and empty at M = (ir, ir). Conse- 
quently, all bands cross the Fermi energy along the TM 
direction, and discontinuinities can be seen at 5 Fermi 
wave vectors k F . Due to symmetries £/ Qr (k F ) may van- 
ish in some cases. For example, at the Fermi wave vector 
kF = (7r/2,7r/2) we have £|| (kF ) = 0, and one can easily 
prove (1, 0, —1, 0, 1) to be an eigenvector of the hopping 
matrix e Q ^(kF) with eigenvalue r\ = 0. This implies 
that there is no discontinuity of n a= 2 (and n a= i) at 
k=(7r/2,7r/2). 

For the interacting system, <5rt Q (k F ) is reduced by 
the layer-dependent factor z a < 1. The Fermi sur- 
faces themselves, however, remain unchanged since 
det(e Q( a(k)) = implies det(y / i^e Q/ 3(k)y / z^') = for 
any k. Generally speaking, the invariance of the Fermi 
surfaces is a consequence of the manifest particle-hole 
symmetry at half-filling and of the local approximation 
for the self-energy. 

From Eq. ([l3]) we have: 

1 f°° 1 

Z » = WT. — Im °t h) (k, E + z0+) dE , (18) 

II k J — oo " 

which shows that z a also yields the weight of the coher- 
ent quasi-particle peak in the layer-resolved density of 
states. z a can serve as an "order parameter" for the 
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Fig. 5. Layer-dependent quasi-particle weight z a as a func- 
tion of U for a d — 5 layer film in sc(100), sc(110) and sc(lll) 
geometry. Calculation for n s = 8 (solid lines) and n s — 10 
(dots). The insets show z a (U) in the critical regime. 



Mott transition. For z a — the system is a Mott- 
Hubbard insulator. 

Fig. 5 shows z a as a function of U in the metallic so- 
lution for the sc(100), sc(110) and the sc(lll) films with 
thickness d = 5. The layer-dependent quasi-particle 
weight decreases from its non-interacting value z a = 1 to 
z a = 0. In the weak-coupling regime there is a quadratic 
dependence 1 — z a (U) oc U 2 which is consistent with per- 
turbation theory in U p6[ . An almost linear dependence 
is seen for the sc(lll) film for a = 1 and a = 3 which 
indicates an early breakdown of perturbation theory in 
this case. For each film geometry considered, there is a 
unique critical interaction U C 2 where all functions z a (U) 
simultaneously approach zero; the whole system can be 
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Fig. 6. Layer dependence of the quasi-particle weight for 
U = 6 and U = 12 in the d = 17 sc(100) film. 



either metallic or insulating. 

Depending on the film geometry, we notice a weak 
(sc(100)) or a rather strong (sc(lll)) layer dependence 
of the quasi-particle weight. For weak and intermediate 
couplings we may observe an oscillatory dependence on 
a (sc(100)). For U closer to U C 2, the behavior changes 
qualitatively; here z a monotonously increases with in- 
creasing distance from a film surface. This is a typical 
effect which is also observed for still thicker films. Fig. 6 
illustrates this fact. It shows a film profile of the quasi- 
particle weight for a d = 17 sc(100) film at U = 6 (z a 
oscillating) and U =12 (z a monotonous). 

In all cases the quasi-particle weight of the surface 
layer z a =\ is significantly reduced compared with a = 
2,3. As has been mentioned above, this can be at- 
tributed to the lowered film-surface coordination num- 
ber. For d = 5 (Fig. 5) the critical interaction U C 2 is 
significantly reduced compared with the bulk (d i— » oo) 
value. For the sc(100) film we find U C 2 = 15.8, for the 
sc(110) U C 2 = 15.4, and for the sc(lll) geometry we 
have U C 2 = 15.2. This has to be compared with the 
bulk value U C 2 = 16.0 (see Fig. 1). The critical interac- 
tion is the lower the more open is the film surface. We 



have n 



(100) 



5, n 



(110) 



4, n 



(in) 



3 for the coordina- 



tion numbers at the film surfaces while in the film center 
n = 6. (Errors for the U C 2 values are discussed below, 
the observed trends are not affected). 




Fig. 7. Layer-dependent 1/E (low-energy) coefficient as a 
function of U for the insulating phase and different film ge- 
ometries. Film thickness: d = 5. 



The same trend is also found for the critical interac- 
tion U c i where the insulating solution disappears. Fig. 7 
shows the U dependence of the 1/E coefficient in the 
low-energy expansion of the self-energy. For the sc(100) 
film we find U cl = 11.3 while U cl = 11.1 for the sc(110) 
geometry. With U c \ = 10.9 the lowest value is ob- 
served for the sc(lll) film. Contrary to the quasi- 
particle weight in the metallic solution, there is always a 
monotonous increase of the 1/E coefficient when passing 
from the film center to one of the surfaces for the insu- 
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Fig. 8. Layer-dependent 1/E coefficient in the low- and the 
high-energy expansion of the self-energy as a function of U. 
Insulator on the d = 5 sc(100) film. 



9 



0.2 



y 1 


1 i 1 i 

metal 

insulator - 


a=1\ 


\\2,3 


■ sc(100) 
d=5 


2.3^v 
i ^ i 



u 

Fig. 9. Layer-dependent double occupancy (n,cr^i-cr) (with 
i 6 a) as a function of U. Metallic (solid lines) and the 
insulating solution (dashed lines) for the d = 5 sc(fOO) film. 



lating solution. In particular, the surface-layer value is 
significantly enhanced. As U i— > U c \ , the 1 /E coefficient 
non-continuously drops to zero. For the thin-film ge- 
ometries this is much more apparent than for the bulk 
(see Fig. 1). 

The almost linear U dependence of the 1/E coefficient 
well above U c i (Fig. 7) alters for still higher U : The 
regime of very strong Coulomb interaction is shown in 
Fig. 8 for the sc(f 00) film. Eventually, for U i— > oo, there 
is a quadratic dependence of the 1/E coefficient on U. 
It smoothly approaches the quadratic U dependence of 
the high-energy expansion of the self-energy: E Q (E) = 

U/2 + U 2 /AE H Q. The system behaves more and 

more atomic-like. 

The layer-dependent average double occupancy is 
shown in Fig. 9 as a function of U. For small U and 
all layers, (n|nj ) decreases linearly as in the Brinkman- 
Rice solution ||17). At U = U C 2 it smoothly joins with 
the (non-zero) average double occupancy of the insu- 
lating solution. For the insulator (n-fnj.) increases with 
decreasing U down to U c \. For both, the metal and 
the insulator, double occupancies are suppressed more 
strongly at the film surfaces compared with the film cen- 
ter. Again, this is due to the stronger effective Coulomb 
interaction U /\fK which results from the reduced vari- 
ance A of the surface free local density of states. 

Within the ED method a second-order transition from 
the metallic to the insulating phase is predicted at U = 
U C 2- This can be seen from the U dependence of the 
respective internal energies: A layer-dependent kinetic 



1.0 




u 

Fig. fO. Kinetic, potential and total energy per site as func- 
tions of U in the metallic (solid lines) and the insulating 
solution (dashed lines) for the d — 5 sc(f 00) film. 

and potential energy per site can be defined as: 

V a = U{n n nn) (i£a). (19) 

The double occupancy (n^n^) and thus the potential 
energy can be obtained from the (ED) solution of the 
impurity models directly. The kinetic energy may be 
calculated from: 

T a = lim ^ReY / {iE n G(iE n )-l) + ^-2U(n i1 n il ). 

(20) 

Fig. 10 shows the kinetic (E kin = {1/ d)Y^, a T a ), po- 
tential (Spot = (V^)Ea^) an d t°t a l energy per site 
(E = E kin + E pot ) for the d = 5 sc(100) film. Due to the 
presence of the film surface, the absolute value of the ki- 
netic energy per site at U = is somewhat lowered com- 
pared with kinetic energy per site of the D = 3 lattice 
which amounts to Elf n ~ 3 ^ = —2.005. Comparing the dif- 
ferent energies in the coexistence region U c \ < U < U C 2, 
we notice that the gain in kinetic energy is almost per- 
fectly outweighted by the loss in potential energy when 
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Fig. 11. U dependence of the quasi-particle weight 23 and 
of the coefficient Q3 in the central layer of the d — 5 sc(100) 
film. Solid lines: n s = 8 (metal) and n s = 7 (insulator). 
Circles: n s = 10 (metal) and n s — 9 (insulator). 



passing from the insulating to the metallic solution at 
a given U. A similar cancellation effect has been ob- 
served beforehand for the D = 00 Bethe lattice fl33| . 
Fig. 10 shows a remaining very small difference in the 
total energy which favors the metallic solution in the 
entire coexistence region. 

Although this is a consistent result within the ED 
method, it must be considered with some care: For the 
metallic solution close to U C 2 the relevant energy scale 
is determined by the width of the quasi-particle reso- 
nance near E = 0. This width is approximately given 
by z a W where W is the free band width. Since z a i— > 
for U 1— ► U C 2, one should expect that finite-size effects 
become important here and that an accurate determina- 
tion of the internal energy cannot be achieved by the ED 
method. The same holds for the determination of U C 2 
(and also of critical exponents). With the ED method 
we cannot access the very critical regime. 

Fig. 11 gives an example. Here we compare the quasi- 
particle weight Z3 at the center of the d — 5 sc(100) film 
obtained for n s — 8 with the result for n s — 10. We no- 
tice that finite-size effects are unimportant for z s > 0.01. 
Generally, n s = 8 sites are sufficient for convergence if U 
is well below U C 2- This can also be seen in Fig. 5 where 
a comparison with the results for n s = 10 is shown at 
a few points. For z < 0.01, however, there are non- 
negligible differences: U C 2 is lowered by about 5% when 
passing from n s — 8 to n s — 10 (Fig. 11). This gives 




Fig. 12. Quasi-particle weight of the central layer a = d/2 
(a — (d + l)/2, respectively) as a function of U for sc(100), 
sc(110) and sc(lll) films (solid lines) with a thickness d 
ranging from d = 1 (two-dimensional Hubbard model) up 
to d = 5 (d = 6, d = 8, respectively) and the result for 
d 1— > 00 (bulk, dashed line) . Note the different U scales. 



an estimate for the error in the determination of U C 2- 
The smallest reliable value for the quasi-particle weight 
Zmin determines the "energy resolution" AE that can 
be achieved by the ED method for a given n s . We have 
AE w z min W. For n s — 8 and with z min ps 0.01 this 
yields AE » 0.1. The error that is introduced by the fi- 
nite low-energy cutoff is of minor importance compared 
with the error due to finite-size effects. This is plausible 
since the smallest fictitious Matsubara frequency can be 
made as small as AE (here we have 7T/3 -1 w 0.05). 

In the particle-hole symmetric case there is always 
one conduction-band energy with = (per layer). 
The corresponding hybridization strength vjf^ vanishes 
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Fig. 13. Thickness dependence of the critical interaction U C 2 
for sc(100), sc(110) and sc(lll) films. 



Fig. 14. Thickness dependence of XJ C \. 



as U i ► U C 2 in the self-consistent calculation. There- 
fore, for the insulator at U > U C 2 we are effectively left 
with n s — 1 sites that are coupled via the hybridization 
term. The metallic solution for n s = 8 (n s = 10) thus 
merges with the insulating solution for n s = 7 (n s = 9) 
at U c2 - Since the quasi-particle resonance is absent in 
the insulating solution, there are no problems with the 
n s convergence in this case. The comparison between 
the results for n s — 7 and n s — 9 is also shown in 
Fig. 11. The critical interaction U c i can be determined 
much more precisely. 

A basic question in a study of the Mott transition 
in thin films concerns the thickness dependence of the 
critical interactions U C 2 and U c \. Since it is the general 
trend that is of primary interest, the mentioned am- 
biguity in the determination of [/ c2 plays a minor role 
only. For all geometries and thicknesses considered, we 
found U C 2 and U c \ to be unique, i. e. for a given sys- 
tem all layer-dependent quasi-particle weights and also 
all 1/E coefficients vanish at a common critical value 
of the interaction, respectively. It is thus sufficient to 
concentrate on the quantities at the film center. 

Fig. 12 shows the central-layer quasi-particle weight z c 
as a function of U for the different systems in the region 
close to U C 2 (d) . It is worth mentioning that for all three 
film geometries and for all d > 2 the functions z c (U) are 
more or less only shifted rigidly to lower interactions 
compared with the result z(U) for the bulk (which of 



course is independent from the choice of the film sur- 
face). This is in clear contrast to the trend of z c (U) for 
weaker interactions (see the result for the sc(lll) film 
in Fig. 5, for example): In the critical regime but also 
at weaker interactions where finite-size effects arc unim- 
portant, the U dependence of z a is rather insensitive to 
details of the film geometry. On the other hand, the 
shift of z c (U) with respect to z{U) docs depend on the 
thickness d. In all cases z c (U) shifts to higher interac- 
tion strengths with increasing d and converges to the 
bulk curve z(U) finally. This also implies that without 
any exception z c increases with increasing d at a given 
U. We can state that the general trends are remarkably 
simple. 

From the zero of z c (U) we can determine the criti- 
cal interaction U C 2 for a given geometry and thickness. 
The results for the different systems are summarized in 
Fig. 13. U c2 is a monotonously increasing function of 
d for the sc(100), the sc(110) as well as for the sc(lll) 
films. For a finite thickness, U C 2(d) is always smaller 
than the bulk value which is apparently approached only 
in the limit d i— > oo. For a given film thickness, U C 2 in- 
creases when passing from the sc(lll) to the sc(110) 
and to the sc(100) film. (Note that for the d = 1 sc(lll) 
film U c2 = since there is no intra-layer hopping in 
this case.) The convergence with respect to d is the 
fastest for the sc(100) and the slowest for the sc(lll) 
films. These trends seem to be related to the differences 
in the film-surface coordination numbers: ri( 100 ) = 5, 
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The interaction strength at which the central-layer 
1/E coefficient a c drops to zero determines C7 cl (see 
Fig. 7). Its thickness dependence is shown in Fig. 14 for 
the different geometries. The result is rather surprising: 
Apart from small uncertainties in the determination of 
U C 2 and U c i , the relative thickness and geometry depen- 
dence of Ud is the same as that of E7 c2 . Only the abso- 
lute values for U c \ are smaller: U c i{d) converges to the 
bulk value U c i = 11.5 which is well below U C 2 = 16.0. 
If we rescale U c i(d) for the different geometries by the 
same (bulk) factor r = U c2 /U c i ~ 1-39, we end up with 
the results for U C 2(d) shown in Fig. 13 within a tolerance 
(|r • U c i(d) - U c2 (d)\/U c2 < 0.005) that is much smaller 
than e. g. the error in the determination of [/ c2 due to 
finite-size effects. The latter is irrelevant here if one 
assumes U C 2 to be overestimated by the same constant 
factor, for the bulk and for the films, which would only 
change the ratio r. Looking at the trends in the results 
for n s = 8 in Fig. 12, this assumption is quite plausible. 
The found relation between U c2 and £7 cl is surprising 
because the disappearance of the insulating solution for 
77 h- » U c i is expected to be of different nature compared 
with the breakdown of the Fermi-liquid metallic phase 
as 77 i ^ U C 2- 

For U C 2, all the details of its geometry and thickness 
dependence can be understood by a simple but instruc- 
tive argument. The main idea has first been developed 
by Bulla (SCj for the D = oo case. Here we discuss a 
generalization for the film geometries: 

The argument assumes that for U > U C 2 one can dis- 
regard the effects of the high-energy charge excitation 
peaks at E » ±77/2 and that the quasi-particle reso- 
nance near E — results from a SIAM hybridization 
function A^(E) = £ fc (V fe (a) )7(.E ~ 4" } ) that consists 
of a single pole at E = 0: 



A^{E) 



\(a) 

E 



(21) 



With the index N we refer to the N-th step in the itera- 
tive solution of the DMFT self-consistency equation. For 
each layer a = 1, d the one-pole structure of A^ a '(E) 
corresponds to an n s — 2 site SIAM with a conduction- 
band energy at = £7/2 and hybridization strength 
V^( Q ) = (A^) 1 / 2 . The analytic solution up to second 
order in /U (cf. Ref. pl[ ) yields two peaks in the 
impurity spectral function at E ss ±[7/2 as well as two 
peaks near E = which build up the Kondo resonance 
for 77 i — ► 77 c2 • The weight of the resonance is thus given 



18(V^) 2 _ 36 (a) 



(22) 



In the self-consistent solution this is also the layer- 
dependent quasi-particle weight which determines the 
coherent part of the film Green function in the low- 
energy regime via Eq. (|l3|). The coherent part of the 
on-site Green function in the a-th layer may be written 
as: 

Gi coh -\E) = z a ■ G a (E) = *" . (23) 

E-M£>F a {E) 

Here, G a (E) is the on-site element of the free film 
Green function, but calculated for the renormalized hop- 
ping matrix e Q/3 (k) i-> ^z^e ct/3 (k) v /z5' (see Eq. (Q). 
The second expression represents the first step in a 
continued-fraction expansion which involves the second 
moment Afl 2) = / dEE 2 G a {E) of G a (E). For the re- 
mainder we have F a (E) =1/E + 0(E- 2 ). 

Starting from Eq. (§|) in the N-th step, the DMFT 
self-consistency equation (Q) yields a new hybridization 
function via A^(E) = E-(e d -fi)-'E a (E)-G a (E)- 1 . 
Using e d = 0, n = U/2, Z a (E) = U/2 + (1 - z^E + ■■■ 
and Eq. (p3|), we get: 

A<> a \E) = —M^F a {E) (24) 

z a 

for energies close to E = 0. Insisting on the one-pole 
structure of A (a \E) = A$ +1 /E for C7 i-> C7 c2 , we must 
have F a (E) = 1/E. This amounts to replacing the co- 
herent part of the film Green function by the simplest 
Green function with the same moments up to the sec- 
ond one. To express in terms of A^ , we still 



"{2) 

need Ma ■ Let us introduce the intra- and inter-layer 
coordination numbers q and p (e.g. q = 4 and p = 1 for 
the sc(100) geometry). The second moment is easily cal- 
culated by evaluating an (anti-)commutator of the form 
([[[c,H ]-,ffo]-,ct] + >. This yields: 

7Wl 2) = z a (qz a + pz a -i + pz a+ i) t 2 . (25) 

We also define the following tridiagonal matrix with di- 
mension d: 

( q P \ 

P 



K([7) = 



36 1 2 
U 2 



V 



P 

q 
p 



(26) 



/ 



Inserting ( |25| ) and ( |22| ) into (^ 
self-consistency equation for 77 





yields a "linearized" 
U c2 : 



N 



(27) 
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Fig. 15. Critical interactions Uc2 an d U^i for sc(100) 
(squares), sc(110) (diamonds) and sc(lll) films (circles) as in 
Fig. 13 and Fig. 14, but plotted against q + 2pcos(7r/(d+l)). 
The dashed lines are linear fits to the data. 



A fixed point of K([7) corresponds to a self-consistent 
solution. Let X r (U) denote the eigenvalues of K(t/). We 
can distinguish between two cases: If |A r (£/)| < 1 for all 
r = 1, d, there is the trivial solution limjvwoo = 
only. This situation corresponds to the insulating so- 
lution for U > U C 2- Contrary, if there is at least one 



X r (U) > 1, diverges exponentially as N i— ► oo. 

This indicates the breakdown of the one-pole model for 
the hybridization function in the metallic solution for 
U < U C 2- The critical interaction is thus determined by 
the maximum eigenvalue: 



A m ax(C^c2) 



1 



(28) 



The eigenvalues of a tridiagonal matrix ( |2q ) can be cal- 
culated analytically for arbitrary matrix dimension d 
]52j : The X r (U) are the zeros of the d-th degree Cheby- 
shev polynomial of the second kind. Solving (|2^) for C/ c2 
then yields: 



U c2 (d,q,p) 




2pcos 



1 



(29) 



Fig. 15 shows that the numerical results for U C 2 can 
be well described by this simple formula. Plotting U 2 2 



against g + 2pcos(7r/(<i-|- 1)), yields a linear dependence 
as a good approximation. The slope in the linear fit 
to the data, however, turns out to be (6.61i) 2 > (Qt) 2 . 
This shows that there is a systematic overestimation of 
the critical interaction by about 10% compared with Eq. 
(29). Partially, this has to be ascribed to finite-size ef- 
fects in the ED method since the error estimated from 
the results in Fig. 11 is of the same order of magnitude. 
(As mentioned before, the fact that the error is a sys- 
tematic one, is explained by the very simple trends seen 
in Fig. 12). On the other hand, we also have to bear in 
mind that ( p9[ ) rests on some simplifying assumptions. 

While we have achieved a satisfactory understanding 
of the results for U C 2, it remains unclear to us why Eq. 
( p9| ) (with the factor 6t replaced by a suitable constant) 
also well describes the thickness and geometry depen- 
dence of U c i. Again, Fig. 15 shows a linear trend. The 
slope is (4.74i) 2 . The above-developed argument, how- 
ever, obviously breaks down. 

It is instructive to compare these results with the anal- 
ogous results for a thermodynamic phase transition in 
thin films. Within the mean-field approach, the layer- 
dependent magnetization m a in Ising films with cou- 
pling constant J is determined by the self-consistency 
equation: 



m 



tanh 



J 



2k B T 



(qm a + pm c 



(30) 



The equation can be linearized for temperatures T near 
the Curie temperature Tc where m a <C 1. Comparing 
with Eq. (|27]) yields the following analogies: 

m a z a , J/2 36i 2 , k B T U 2 , k B T c & U 2 2 . 

(31) 

The exact Curie temperature in thick (e? i— > oo) Ising 
films is expected |54| to obey the power law (T c (oo) — 
T c (d))/T c (p6) — CocP^ where the shift exponent A = 
\jv is related to the (D = 3) critical exponent v of the 
correlation function. Within the mean-field approach to 
the magnetic phase transition in Ising films one obtains 
A = 2. 

The same exponent is found within dynamical mean- 
field theory applied to the Mott transition in Hubbard 
films. Expanding Eq. (|29[) in powers of 1/d, we obtain: 



U c2 - U c2 (d,q,p) 



U, 



<-2 



C (q,p) ■ d 



-A 



(32) 



where U C 2 denotes the bulk value, Co(q,p) 
71-2 q/(q + 2p) and A = 2 the exponent. 
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IV. CONCLUSION 

By the invention of dynamical mean-field theory we 
are in a position to treat itinerant-electron models on 
the same footing as does the Weiss molecular-field the- 
ory with respect to localized spin models. To gain a first 
insight into collective magnetism for thin-film and sur- 
face geometries, the molecular-field approach has suc- 
cessfully been applied to the Ising or the Heisenbcrg 
model in the past. It is well known that the reduced 
system symmetry may result in characteristic modifica- 
tions of the magnetic properties which are interesting on 
their own but also with respect to technical applications. 
The Weiss theory is able to describe a major number of 
magnetic properties qualitatively correct. 

To study the correlation-induced Mott transition from 
a paramagnetic metal to a paramagnetic insulator, we 
need to consider an itinerant-electron model. Presum- 
ably, the Hubbard model is the simplest one in this re- 
spect. Just as the Weiss theory for the magnetic prop- 
erties of Ising and Heisenberg films, the DMFT provides 
the first step in an understanding of the Mott transition 
and the related electronic properties in Hubbard films. 
This has been the central idea of the present study. We 
have generalized the mean-field equations for the appli- 
cation to systems with reduced translational symmetry 
and have solved them using the exact-diagonalization 
method. Let us briefly list up the main results found: 

Similar to the results for the infinite-dimensional 
Bethe lattice, we find a metallic and an insulating so- 
lution of the mean-field equations at T = for all film 
geometries considered. These coexist in a certain range 
of the interaction strength. The metallic solution is sta- 
ble against the insulating one in the whole coexistence 
region (however, in the critical regime the ED approach 
is questionable). 

Generally, the breakdown of translational symmetry 
with respect to the film normal direction leads to mod- 
ifications of the electronic structure. For film thickness 
d i ► oo these disappear in the bulk, while surface effects 
are still present. The finite film thickness as well as the 
presence of the surface manifest themselves in a signifi- 
cant layer dependence of the on-site Green function and 
the self-energy for both, the metallic and the insulating 
solution. The layer dependence is found to be the more 
pronounced the more open is the film surface. In thicker 
films surface effects quickly diminish when passing from 
the top layer to the film center. 

In particular, we have considered the so-called layer- 
dependent quasi-particle weight z a for the metallic 
phase as a function of U, d and geometry. For U 7^ 0, 



z a < 1 is the reduction factor of the discontinuous drop 
of the momentum-distribution function in the a-th layer 
at each of the d one-dimensional Fermi- "surfaces" or, 
equivalently, the weight of the coherent quasi-particle 
peak in the local density of states of the a-th layer. 
In all cases the quasi-particle weight at the film sur- 
faces is found to be significantly reduced which is due 
to the surface enhancement of the effective Coulomb 
interaction U / VA (A is the variance of the free local 
density of states). At the film surfaces the electrons 
are "heavier", double occupancies are suppressed more 
effectively The layer dependence of the quasi-particle 
weight is oscillatory in some cases at small and interme- 
diate interaction strengths. As U approaches the critical 
regime, the behavior changes qualitatively; here z a is al- 
ways monotonously increasing with increasing distance 
from a film surface. Furthermore, for fixed U, z a always 
increases with increasing d - the general trends are re- 
markably simple. 

The low-energy electronic structure in the insulating 
solution is governed by the l/E coefficient in the low- 
energy expansion of the self-energy. Generally, the layer 
dependence of the coefficient is less spectacular com- 
pared with the layer dependence of z a in the metallic 
solution. At the film surfaces the 1 /E coefficient is some- 
what enhanced and in all cases monotonously decreases 
with increasing distance from a surface. 

For a given film geometry and thickness there is 
a unique critical interaction C/ c2 at which all (layer- 
dependent) effective masses z^ 1 diverge. This implies 
that the whole film is either metallic or a Mott insula- 
tor. For all cases investigated, we did not find a surface 
phase that differs from the bulk phase. This may change 
if non-uniform model parameters, e. g. a modified sur- 
face U, are considered. The question is left for future 
investigations. 

While a precise determination of the critical inter- 
action is not possible by means of the ED approach, 
general trends can be derived safely. It is found that 
U C 2 is strongly geometry and thickness dependent. It 
monotonously increases with increasing film thickness 
and smoothly approaches the bulk value from below. 
For finite d, U C 2 is the smaller and the convergence to 
the bulk value is the slower the more open is the film 
surface. The same trends are also seen for the criti- 
cal interaction U& at which the insulating solution dis- 
appears. In fact, within numerical accuracy a simple 
relation, U C 2{d) = r ■ U c i(d), seems to hold. The geome- 
try and thickness dependence of U C 2 can be understood 
qualitatively by (approximately) linearizing the mean- 
field equations for U 1— > After proper rescaling, 
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the resulting analytical expression for U C 2 (d, q, p) can 
even quantitatively reproduce the numerical data. The 
argument, however, is not applicable to the critical in- 
teraction U c i and thus cannot explain the "empirically" 
found relation between f7 cl and U C 2- An effective, lin- 
earized theory that is valid for U ► U c \ remains to be 
constructed. Calculating the Curie temperature of Ising 
films within the molecular-field theory, yields exactly the 
same trends. This analogy suggests that the observed 
geometry and thickness dependence may be typical for 
the mean-field treatment of the problem. Whether or 
not there are qualitative changes if U c i^{d,q,p) could 
be calculated beyond mean-field theory, remains to be 
another open problem. 

Finally, one could think of an experimental investiga- 
tion of the Mott transition in thin films. The present 
study shows that the strong thickness dependence of 
the critical interaction has its origin in the reduced 
coordination number at the film surface. A monotonous 
increase of the critical interaction with increasing d 
should thus be expected also for temperatures above 
the Neel temperature where we have a magnetically 
disordered state. Studying the Mott transition in a 
bulk material, one needs to "vary U" (or temperature) 
experimentally. Contrary, for a thin-film geometry the 
transition may take place also by varying d. If it is 
possible to grow crystalline films of a metallic material 
which in the bulk is close to the Mott transition, one 
may observe insulating behavior in ultrathin films and a 
transition to a metallic phase with increasing thickness. 
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